/*
    Copyright (C) 2009, 2010, 2011 Daniel Gutson and Martin Gioiosa, FuDePAN

    This file is part of Biopp.

    Biopp is free software: you can redistribute it and/or modify
    it under the terms of the GNU General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.

    Biopp is distributed in the hope that it will be useful,
    but WITHOUT ANY WARRANTY; without even the implied warranty of
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
    GNU General Public License for more details.

    You should have received a copy of the GNU General Public License
    along with Biopp.  If not, see <http://www.gnu.org/licenses/>.

    NOTE: This file is in prototype stage, and is under active development.
*/

#include <iostream>
#include "biopp/biopp.h"
#include "biopp/bio_molecular/biopp_debug.h"

using namespace std;

namespace biopp
{

GeneticCode GeneticCode::instance;

GeneticCode::GeneticCode()
{
    add("ttt", "F"),   add("tct", "S"), add("tat", "Y"),   add("tgt", "C"),
        add("ttc", "F"),   add("tcc", "S"), add("tac", "Y"),   add("tgc", "C"),
        add("tta", "L"),   add("tca", "S"), add("taa", Aminoacid::STOP_CODON), add("tga", Aminoacid::STOP_CODON),
        add("ttg", "L"),   add("tcg", "S"), add("tag", Aminoacid::STOP_CODON), add("tgg", "W"),
        add("ctt", "L"),   add("cct", "P"), add("cat", "H"),   add("cgt", "R"),
        add("ctc", "L"),   add("ccc", "P"), add("cac", "H"),   add("cgc", "R"),
        add("cta", "L"),   add("cca", "P"), add("caa", "Q"),   add("cga", "R"),
        add("ctg", "L"),   add("ccg", "P"), add("cag", "Q"),   add("cgg", "R"),
        add("att", "I"),   add("act", "T"), add("aat", "N"),   add("agt", "S"),
        add("atc", "I"),   add("acc", "T"), add("aac", "N"),   add("agc", "S"),
        add("ata", "I"),   add("aca", "T"), add("aaa", "K"),   add("aga", "R"),
        add("atg", "M"),   add("acg", "T"), add("aag", "K"),   add("agg", "R"),
        add("gtt", "V"),   add("gct", "A"), add("gat", "D"),   add("ggt", "G"),
        add("gtc", "V"),   add("gcc", "A"), add("gac", "D"),   add("ggc", "G"),
        add("gta", "V"),   add("gca", "A"), add("gaa", "E"),   add("gga", "G"),
        add("gtg", "V"),   add("gcg", "A"), add("gag", "E"),   add("ggg", "G")
        ;
}

ExtendedGeneticCode ExtendedGeneticCode::instance;

ExtendedGeneticCode::ExtendedGeneticCode()
{
    a[0] = PseudoNucleotide::A;
    a[1] = PseudoNucleotide::STOP_PN;
    t[0] = PseudoNucleotide::U;
    t[1] = PseudoNucleotide::STOP_PN;
    c[0] = PseudoNucleotide::C;
    c[1] = PseudoNucleotide::STOP_PN;
    g[0] = PseudoNucleotide::G;
    g[1] = PseudoNucleotide::STOP_PN;

    r[0] = PseudoNucleotide::A;
    r[1] = PseudoNucleotide::G;
    r[2] = PseudoNucleotide::STOP_PN;
    y[0] = PseudoNucleotide::U;
    y[1] = PseudoNucleotide::C;
    y[2] = PseudoNucleotide::STOP_PN;
    k[0] = PseudoNucleotide::U;
    k[1] = PseudoNucleotide::G;
    k[2] = PseudoNucleotide::STOP_PN;
    m[0] = PseudoNucleotide::A;
    m[1] = PseudoNucleotide::C;
    m[2] = PseudoNucleotide::STOP_PN;
    s[0] = PseudoNucleotide::C;
    s[1] = PseudoNucleotide::G;
    s[2] = PseudoNucleotide::STOP_PN;
    w[0] = PseudoNucleotide::A;
    w[1] = PseudoNucleotide::U;
    w[2] = PseudoNucleotide::STOP_PN;

    b[0] = PseudoNucleotide::U;
    b[1] = PseudoNucleotide::C;
    b[2] = PseudoNucleotide::G;
    b[3] = PseudoNucleotide::STOP_PN;
    d[0] = PseudoNucleotide::A;
    d[1] = PseudoNucleotide::U;
    d[2] = PseudoNucleotide::G;
    d[3] = PseudoNucleotide::STOP_PN;
    h[0] = PseudoNucleotide::A;
    h[1] = PseudoNucleotide::U;
    h[2] = PseudoNucleotide::C;
    h[3] = PseudoNucleotide::STOP_PN;
    v[0] = PseudoNucleotide::A;
    v[1] = PseudoNucleotide::C;
    v[2] = PseudoNucleotide::G;
    v[3] = PseudoNucleotide::STOP_PN;

    n[0] = PseudoNucleotide::A;
    n[1] = PseudoNucleotide::U;
    n[2] = PseudoNucleotide::C;
    n[3] = PseudoNucleotide::G;
    n[4] = PseudoNucleotide::STOP_PN;

    expansions[0] = &a[0];
    expansions[1] = &t[0];
    expansions[2] = &c[0];
    expansions[3] = &g[0];
    expansions[4] = &r[0];
    expansions[5] = &y[0];
    expansions[6] = &k[0];
    expansions[7] = &m[0];
    expansions[8] = &s[0];
    expansions[9] = &w[0];
    expansions[10] = &b[0];
    expansions[11] = &d[0];
    expansions[12] = &h[0];
    expansions[13] = &v[0];
    expansions[14] = &n[0];
}

TripletsIterator GeneticCode::triplets(Aminoacid amino) throw(InvalidAminoacid)
{
    typedef std::pair < std::multimap<Aminoacid, Triplet>::const_iterator,
            std::multimap<Aminoacid, Triplet>::const_iterator > piter;

    class GenCodeTripletsIterator :  public TripletsIteratorInterface
    {
        const piter pit;
        std::multimap<Aminoacid, Triplet>::const_iterator it;

        virtual void restart()
        {
            it = pit.first;
        }

        virtual bool end() const
        {
            return it == pit.second;
        }

        virtual void next()
        {
            ++it;
        }

        virtual void endOfUse()
        {
            delete this;
        }

        virtual const Triplet& get() const
        {
            return it->second;
        }

        virtual GenCodeTripletsIterator* clone() const
        {
            return new GenCodeTripletsIterator(pit);
        }

    public:
        GenCodeTripletsIterator(piter pit)
            : pit(pit)
        {
            it = pit.first;
        }
    };

    piter pit = instance.a2t.equal_range(amino);

    if (pit.first == pit.second)
        throw InvalidAminoacid();

    return TripletsIterator(new GenCodeTripletsIterator(pit));
}

void Triplet::from_string(const std::string& s) throw(InvalidNucleotide)
{
    if (s.size() != 3)
        throw InvalidNucleotide("Triplet size invalid:" + mili::to_string(s.size()));

    triplet[0] = s[0];
    triplet[1] = s[1];
    triplet[2] = s[2];
}

TripletsIterator Triplet::getSiblingTriplets() const
{
    class SiblingsTripletsIterator : public TripletsIteratorInterface
    {
        const Triplet orig;
        TripletsIterator it;

        virtual void restart()
        {
            it.restart();

            if (*it == orig)
                ++it;
        }

        virtual bool end() const
        {
            return it.end();
        }

        virtual void next()
        {
            it.next();
            if (*it == orig)
                it.next();
        }

        virtual void endOfUse()
        {
            delete this;
        }

        virtual const Triplet& get() const
        {
            return *it;
        }

        virtual SiblingsTripletsIterator* clone() const
        {
            return new SiblingsTripletsIterator(orig);
        }
    public:
        SiblingsTripletsIterator(const Triplet& orig)
            : orig(orig)
        {
            it = GeneticCode::triplets(GeneticCode::translate(orig));

            if (*it == orig)
                ++it;
        }
    };

    return TripletsIterator(new SiblingsTripletsIterator(*this));
}

size_t Triplet::min_distance(Aminoacid aminoacid) const throw(InvalidAminoacid)
{
    size_t min_dist(4);

    for (TripletsIterator triplet = GeneticCode::triplets(aminoacid); !triplet.end(); ++triplet)
        mili::minimize(min_dist, TripletBase<Nucleotide>::distance(*triplet));

    return min_dist;
}

TripletsIterator GeneticCode::allTriplets()
{
    class AllTripletsIterator : public TripletsIteratorInterface
    {
        std::map<Triplet, Aminoacid>::const_iterator it;

        virtual void restart()
        {
            it = GeneticCode::instance.t2a.begin();
        }

        virtual bool end() const
        {
            return it == GeneticCode::instance.t2a.end();
        }

        virtual void next()
        {
            ++it;
        }

        virtual void endOfUse()
        {
            delete this;
        }

        virtual const Triplet& get() const
        {
            return it->first;
        }

        virtual AllTripletsIterator* clone() const
        {
            return new AllTripletsIterator;
        }

    public:
        AllTripletsIterator()
        {
            restart();
        }
    };

    return TripletsIterator(new AllTripletsIterator);
}

AminoacidsIterator GeneticCode::allAminoacids()
{
    class AllAminoacidsIterator : public AminoacidsIteratorInterface
    {
        std::set<Aminoacid>::const_iterator it;

        virtual void restart()
        {
            it = GeneticCode::instance.aminoacids.begin();
        }

        virtual bool end() const
        {
            return it == GeneticCode::instance.aminoacids.end();
        }

        virtual void next()
        {
            ++it;
        }

        virtual void endOfUse()
        {
            delete this;
        }

        virtual const Aminoacid& get() const
        {
            return *it;
        }

        virtual AllAminoacidsIterator* clone() const
        {
            return new AllAminoacidsIterator;
        }
    public:
        AllAminoacidsIterator()
        {
            restart();
        }
    };

    return AminoacidsIterator(new AllAminoacidsIterator);
}

TripletsIterator Triplet::mutate_single(Nucleotide from, Nucleotide to) const
{
    class MutateSingleIterator : public TripletsIteratorInterface
    {
        const Triplet orig;
        const Nucleotide from;
        const Nucleotide to;

        size_t position;
        mutable Triplet mutated;

        void advance()
        {
            while (position < 3 && orig[position] != from)
                position++;
        }

        virtual void restart()
        {
            position = 0;
            advance();
        }

        virtual bool end() const
        {
            return position > 2;
        }

        virtual void next()
        {
            position++;
            advance();
        }

        virtual void endOfUse()
        {
            delete this;
        }

        virtual const Triplet& get() const
        {
            mutated = orig;
            mutated[position] = to;
            return mutated;
        }

        virtual MutateSingleIterator* clone() const
        {
            return new MutateSingleIterator(orig, from, to);
        }
    public:
        MutateSingleIterator(const Triplet& orig, Nucleotide from, Nucleotide to)
            : orig(orig), from(from), to(to)
        {
            restart();
        }
    };

    return TripletsIterator(new MutateSingleIterator(*this, from, to));
}

void PseudoTriplet::from_string(const std::string& s) throw(InvalidNucleotide)
{
    if (s.size() != 3)
        throw InvalidNucleotide("Triplet size invalid:" + mili::to_string(s.size()));

    triplet[0] = s[0];
    triplet[1] = s[1];
    triplet[2] = s[2];
}

bool PseudoTriplet::isPureTriplet() const
{
    bool isPure = true;
    const PseudoNucleotide last_nucleotide('G');

    size_t i = 0;

    do
    {
        if (triplet.n[i] > last_nucleotide)
            isPure = false;
        ++i;
    }
    while (isPure && i < 3);

    return isPure;
}

TripletsIterator PseudoTriplet::getPureTriplets() const
{
    return ExtendedGeneticCode::possibleTriplets(*this);
}

TripletsIterator ExtendedGeneticCode::possibleTriplets(const PseudoTriplet& pseudo)
{
    typedef std::vector<Triplet> vt;

    class  PseudoTripletsToTripletsIterator :  public TripletsIteratorInterface
    {
        const PseudoTriplet orig;

        vt ctriplets;
        vt::const_iterator it;

        virtual void restart()
        {
            it = ctriplets.begin();
        }

        virtual bool end() const
        {
            return it == ctriplets.end();
        }

        virtual void next()
        {
            ++it;
        }

        virtual void endOfUse()
        {
            delete this;
        }

        virtual const Triplet& get() const
        {
            return *it;
        }

        virtual PseudoTripletsToTripletsIterator* clone() const
        {
            return new PseudoTripletsToTripletsIterator(orig);
        }

        void combinatory(PseudoTriplet& orig, Triplet& aux, size_t n = 0)
        {
            for (const PseudoNucleotide* p = instance.expansions[orig[n].value]; *p != PseudoNucleotide::STOP_PN; ++p)
            {
                aux[n] = *p;

                if (n != 2)
                    combinatory(orig, aux, n + 1);
                else
                    ctriplets.push_back(aux);
            }
        }

        void combinatory(PseudoTriplet& orig)
        {
            Triplet aux;
            combinatory(orig, aux);
        }

    public:
        PseudoTripletsToTripletsIterator(PseudoTriplet pseudo)
            : orig(pseudo)
        {
            combinatory(pseudo);

            it = ctriplets.begin();
        }
    };

    return TripletsIterator(new  PseudoTripletsToTripletsIterator(pseudo));
}
}
